Method for correcting the phase of electromagnetic data

ABSTRACT

Method for identifying, determining and correcting source-related phase errors in data from a controlled source electromagnetic survey by using data from ordinary survey receivers, i.e. without benefit of source monitoring data. Abrupt anomalies indicating source malfunctions are identified ( 71 ) in the time domain by plotting time intervals between neighboring zero crossings or by zero-lag cross correlation between consecutive bins of receiver data, and the amount of the time error ( 73 ) can be determined by performing cross correlation between two bins on either side of an anomaly. In the frequency domain, transmitter anomalies can be identified by looking for discontinuities in plots of phase vs. offset, and the corrective phase shift can be determined by matching the phase on one side of the anomaly to that on the other side. A global time/phase shift ( 76 ) can be determined by using phase frequency-scaling behavior at near offsets.

This application is a National Stage entry under 35 U.S.C. 371 ofPCT/US2007/013737 that published as WO 2008/013609 and was filed on Jun.12, 2007 and claims the benefit of U.S. Provisional application60/833,088 which was filed on Jul. 25, 2006.

FIELD OF THE INVENTION

This invention relates generally to the field of geophysicalprospecting, and more particularly to electromagnetic prospecting forhydrocarbons. Specifically, the invention is a method for detecting,determining, and correcting phase errors in controlled sourceelectromagnetic (“CSEM”) data by using normal-receiver data alone.

BACKGROUND OF THE INVENTION

The CSEM survey technique is an important geophysical tool forprospecting for hydrocarbons in the earth's subsurface. In a CSEMsurvey, an electromagnetic-wave source (transmitter) generates anelectromagnetic wave. The electromagnetic signal induced in the earth bythe transmitter is recorded continuously in time by one or morereceivers. The electromagnetic signal at a receiver location depends onphysical properties, especially the electrical properties, of the mediumthrough which the electromagnetic wave has passed from the source to thereceiver. The behavior of this signal as a function of frequency andtransmitter location or separation (offset) between transmitter andreceiver can be used to estimate the spatially varying resistivity modelof the subsurface. This estimated resistivity model is used foridentifying hydrocarbons in the earth's subsurface.

In a typical marine CSEM survey, a constantly activeelectromagnetic-wave transmitter is towed along a line 11 aboveelectromagnetic receivers 12 deployed on the seafloor 13, as illustratedin FIG. 1. For more details see for example Chapter 12, page 931 inInvestigations In Geophysics No. 3, Electromagnetic Methods In AppliedGeophysics, volume 2, edited by Misac N. Nabighian, Society ofExploration Geophysicists (1991). The receivers typically have multiplesensors designed to record different vector components of the electricand/or magnetic fields. The directions of the received data componentsby the receivers are indicated in FIG. 1 by arrows. Transmitterlocations are illustrated by arrows along line 11 above the receivers.Normally, the separation between neighboring transmitter locations ismuch smaller than that between neighboring receivers.Transmitter-receiver offset is defined as the distance 14 between atransmitter location and a receiver location. The transmitting andreceiving systems typically operate independently (without anyconnection between them), so that receiver data must be synchronizedwith shipboard measurements of the instantaneous transmitter positionand electric current in the transmitter antenna by comparing clock timeson the receivers to time from a shipboard or GPS (Global PositioningSystem) standard.

Every receiver records the electromagnetic signal continuously in timeduring a survey. The data recorded by one sensor (channel) on a receiverare called a common-receiver gather. This gather represents theelectromagnetic signal at the receiver location induced by the source atall different source locations, or at different times during the survey.FIGS. 2A-C depict the measured electrical signals recorded by the threesensors of a receiver within a short time window (only a few transmitterwaveforms) taken out of a common-receiver gather together with (FIG. 2D)the transmitter waveform that gave rise to it. FIG. 2A shows themeasured x-component of the electric field, FIG. 2B the measuredy-component and FIG. 2C the measured z-component. FIG. 2D shows thetransmitter current signal of a square waveform with designed halfperiod 4 s. An anomalous transmitter pulse 21 of 6 s width shows up onall three channels of the receiver. In practice, the transmitter signalmay be a more complex waveform than the square wave depicted in FIG. 2D.

Marine CSEM data are typically interpreted in the temporal frequencydomain. The resulting transformed data will have components at certainfrequencies, determined by the frequency spectrum of the particularsource waveform. After taking out the frequency-dependent effects of thesource and the receiver itself, the signal at a frequency represents theresponse of the earth to an electromagnetic signal at that temporalfrequency. Like any other type of wave, the electromagnetic signal in aCSEM survey has two attributes, amplitude and phase. The signals aretherefore conveniently represented as complex numbers in eitherrectangular (real-imaginary) or polar (amplitude-phase) form.

In practice, the receiver data are usually converted to temporalfrequency by dividing (or “binning”) the recorded time-domain data 31into time intervals, e.g., bins x₁, x₂, and x₃ in FIG. 3A, anddetermining the spectrum within each bin by standard methods based onthe Fourier Transform. FIG. 3B shows the amplitudes of the spectralcomponents from the bin x₃. A typical bin length can be one or severalperiods of the transmitter waveform 32. Some methods of transformingdata to the time-frequency domain include the Short-Time FourierTransform (J. Allen, L. Rabiner, “A Unified Approach to Short-TimeFourier Analysis and Synthesis,” Proc. of the IEEE 65, 1558-64, (1977));the Wavelet Transform (W. C. Lang and K. Forinash, “Time-frequencyanalysis with the continous wavelet transform,” Am. J. Phys. 66,794-797, (1998)); the Wigner-Ville transform (E. Wigner, “On the quantumcorrection for thermodynamic equilibrium,” Phys. Rev. 40, 749-759,(1932)); the Choi-Williams transform (H. Choi and W. Williams, “Improvedtime-frequency representation of multicomponent signals usingexponential kernels,” IEEE Trans. on Acoust., Speech, and SignalProcessing 37, 862-871, (1989)); and the Bessel method (Z. Guo, L. G.Durand, and H. C. Lee, “The time-frequency distributions ofnonstationary signals based on a Bessel kernel,” IEEE Trans. on SignalProc. 42, 1700-1707, (1994)). The present invention is not limited toany particular method or methods for spectral decomposition of CSEM datato the temporal-frequency domain. In the temporal-frequency domain,signals, including both amplitude and phase, of each of thetemporal-frequency components are functions of a bin designationvariable (label).

With each temporal bin is associated a time, typically the Julian dateat the center of the bin. Since the transmitter location is known as afunction of time, these bins may be interchangeably labeled in severaldifferent ways: by Julian date of the bin center; by transmitterposition; by the signed offset distance between source and receiver; bythe cumulative distance traveled by the transmitter relative to somearbitrarily chosen starting point; or sometimes simply by a sequentialbin number.

In many examples of CSEM hardware, data cannot be effectively recordedat the nearest offsets because the dynamic range of the receiver'sdigitizers is too small to accommodate the large dynamic range of thedata. This region is sometimes known as the “saturation zone” andtypically encompasses source-receiver offsets of less than about 500meters. At far offsets, the controlled-source signal is too weak to beobserved above the earth's magnetotellurics energy and other noises suchas those induced by oceanic currents. This signal level may be calledthe “noise floor”.

Both amplitudes and phased of the CSEM data depend on the subsurfaceconductivities. Both of them need to be determined accurately in orderto distinguish characteristics that signal the presence or absence ofconductivity anomalies in the subsurface. Reliable phases of CSEM dataare essential in applying inversion techniques for conductivity-anomalydetections. However, in a CSEM survey, the phase variations of CSEM datacan be caused by many factors that are not related to the variations ofthe subsurface conductivity in the target zones. Some of these factorsare identified in the following list:

(a) the transmitter current waveform deviates from the desired waveform,i.e., transmitter instability;

(b) the transmitter monitoring system fails to measure and/or report theactual transmitter current waveform accurately;

(c) transmitter and receiver signals are recorded separately usingdifferent time bases (clocks) that are not synchronized against a commonGPS time base;

(d) the frequency-dependent responses of a receiver, such as receiveramplifiers, receiver antennae, and their combination, are not calibratedaccurately;

(e) localized changes in the earth's resistivity close to a receivermask the electromagnetic signals from the desired target zones.

Phase errors of the types of categories (a) and (b) above are common,and can occur in both land and marine CSEM surveys. As an example of thetype of problem of category (c) above, if the time origin of thetransmitter waveform (sometimes known as the “initial transmitterphase”) is in error by Δt, the phase of the receiver data in frequencydomain at angular frequency ω will be incremented by an amount ωΔt.Similar factors arise on land, although it is easier to connect both thesource and receivers to a common time reference, thus reducing the phaseerror associated with the clock synchronization.

Some techniques have been disclosed in the published literature forreducing the phase errors caused by some of the factors listed above.Direct measurement of the receiver-clock drift (time error) relative toa time reference (such as GPS) at the start and end of the survey allowsusers to stretch or compress measured data to an estimate of thereference time (Constable et al., “Marine magnetotellurics for petroleumexploration Part 1: A sea-floor equipment system,” Geophysics 63,816-825 (1998)). This correction can reduce the error caused by thereceiver clock drifts. Laboratory measurements of the response of thereceiver's amplifier-antenna system have been used to compensate fieldCSEM data (Ellingsrud, et al., “Remote sensing of hydrocarbon layers byseabed logging (SBL): Results from a cruise offshore Angola,” TheLeading Edge 21, 972-982 (2002)).

As pointed out above, if the time origin of the transmitter waveform isin error by Δt, the phase of the receiver data in frequency domain atangular frequency ω will be incremented by an amount ωΔt. This type oferror can also be caused by undesired variation of the transmitterwaveform. To resolve the problems of transmitter instability and timingerror associated with the transmitter (initial phase error), thetransmitter current is usually recorded continuously during the surveyby the source monitoring system. An independent monitoring receiver hasalso been mounted to the towed underwater transmitter previously tomonitor the transmitter current that is actually injected into the water(MacGregor et al., “The RAMESSES experiment-III. Controlled-sourceelectromagnetic sounding of the Reykjanes Ridge at 57° 45′ N,” Geophys.J. Int. 135, 773-789 (1998)). In all cases, the transmitter signalsrecorded by the monitoring system in the time domain are used toidentify and compensate the phase shifts caused by anomalous transmittersignal variations during the survey.

In principle, a source deconvolution technique can be used to obtain theresponse of the subsurface at each frequency once the source signal isknown reliably even though the source may exhibit undesired variations(Ö. Yilmaz, Seismic Data Processing, Vol. 2 in Investigations inGeophysics, E. B. Neitzel, ed., Society of Geophysics, 498-506 (1987)).However, the deconvolution technique strongly depends on the fidelity ofthe recorded source signals during the survey. This dependence combinedwith the shortcomings of the deconvolution technique itself makes theapplication of this technique to CSEM data difficult. Instead, in CSEMdata processing, one continuous transmitter towline is normally splitinto multiple sublines. Time windows of erratic variation in transmittersignal are skipped. Each subline uses a different transmitter initialtime (or initial phase), determined from the transmitter currentsignals, for processing. The major shortcoming of this method and thedeconvolution technique is the dependence on the measured transmittersignals. The measured transmitter signals may not be reliable and/oravailable due to malfunctions or failures of the source monitoringsystem. Also, splitting one towline into multiple sublines renders CSEMdata processing tedious and error-prone. After data processing with thetransmitter initial times (or phases), normally the data still haveresidual phase errors that need to be corrected by other methods.

Another method to detect and correct phase errors is to determine anytiming error for one common-receiver gather by utilizing the frequencyscaling behavior of the electromagnetic field in a uniform medium (PCTInternational Patent Application No. PCT/US06/46329). The phase of theelectromagnetic field from a dipole in a uniform medium is a function of(R√{square root over (ω)}), where R is the distance between thetransmitter and the receiver (called offset), and ω is the angularfrequency. In a homogeneous-earth model, the phase variations with thefrequency-scaled transmitter-receiver offset, (R√{square root over(ω)}), are the same for different frequencies. In a one-dimensionallayered-earth model, data at small transmitter-receiver offsets shows asimilar scaling behavior. Separation of the phase versus scaled-offsetcurves at different frequencies within a small offset range (e.g. within2 km but outside the saturation zone) indicates a time error assumingthe conductivity variation within the small offset range is small. Atiming error can be determined by maximizing the overlap of the phaseversus scaled-offset curves of different frequencies at near offset. Thetiming error determined from this method includes components from bothtransmitter and receiver clocks. This method is effective in determiningthe global phase shift for a common-receiver gather, ensuring that thephase goes to zero as the source and receiver offset goes to zero, asrequired by the laws of physics. However, this method is not valid whennear-offset data (such as those with offset less than 1 km) do not existor the undesired transmitter variations happen at locations too far fromthe receiver location for the frequency-scaling behavior of the phase tobe valid.

What is needed is a method, applicable to both onshore and offshore CSEMsurveys, that can detect, determine, and correct phase errors related tothe source (i.e., categories (a) and (b) in the list above) based ondata from normal, ordinary survey receivers, i.e. without using a sourcesignal measured by a source monitoring system, so that data phase errorscan be corrected when the source monitoring data are either not reliableor not available. The present invention satisfies this need.

SUMMARY OF THE INVENTION

In one of its embodiments, the present invention is a method forcorrecting source-related phase errors in electromagnetic data signalsrecorded by receivers in a controlled-source electromagnetic survey of asubsurface region, comprising:

(a) selecting a source line from the survey and an ordinary surveyreceiver associated with the selected source line, and obtaining theelectromagnetic data measured by the selected receiver for sourcepositions along the selected source line, said obtained data beingreferred to as a “receiver gather”;

(b) finding a source-related anomaly in the receiver gather byexamination of said data;

(c) determining from the receiver gather a time shift or correspondingphase shift due to the anomaly; and

(d) using the time/phase shift to adjust the phase of the data in thereceiver gather.

In typical applications, the above-described process is repeated forother source-related anomalies in the data record for the selectedreceiver and source line, then repeated for other receivers on thesource line and finally all the preceding is repeated for other sourcelines in the survey. Moreover, some embodiments of the invention includedetermining whether the receiver gather has a global time/phase error,and if so, determining the error and correcting the gather by applying atime/phase shift.

BRIEF DESCRIPTION OF THE DRAWINGS

The present invention and its advantages will be better understood byreferring to the following detailed description and the attacheddrawings in which:

FIG. 1 is a schematic diagram of a marine CSEM purvey withmulti-component receivers located at the sea floor;

FIGS. 2A-D are plots showing the signals recorded by the three electricfield channels (E_(X) E_(Y), and E_(Z)) and the square wave transmittersignal that produced the receiver responses;

FIGS. 3A-B illustrate the process of binning a receiver signal in timeand determining the frequency spectrum within each bin by Fourieranalysis;

FIGS. 4A-C show plots of time intervals between neighboringzero-crossings determined from receiver data for a transmitter ofsquare-wave waveform of half period 4 s as shown in FIG. 2D;

FIGS. 5A-D illustrate the detection of transmitter anomalies anddetermination of time shifts across anomalies using cross correlation;

FIGS. 6A-B are schematic plots showing the time (or phase) shiftdetermination by matching the phase at one side of a transmitter-signalanomaly to that at another side; and

FIG. 7 is a flow chart showing basic steps in one embodiment of thepresent inventive method.

The invention will be described in connection with its preferredembodiments. However, to the extent that the following description isspecific to a particular embodiment or a particular use of theinvention, this is intended to be illustrative only, and is not to beconstrued as limiting the scope of the invention. On the contrary, it isintended to cover all alternatives, modifications and equivalents thatmay be included within the spirit and scope of the invention, as definedby the appended claims.

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

The present invention is a method for identifying, determining, andcorrecting phase errors of CSEM data resulting from transmitteranomalies by using ordinary-receiver data alone. The ordinary receiversin this invention are those receivers deployed for recording the inducedelectromagnetic field by the active source in a CSEM survey (generallyreferred to as the survey receivers, or simply receivers), in contrastto a receiver associated with the source-signal monitoring systemmeasuring the electric current of the source or the directelectromagnetic field from the source. The receiver in thesource-monitoring system is usually connected to the source so that theoffset does not vary during the survey. One significant aspect of thisinvention is its capability to identify, determine, and correct CSEMdata phase errors by using data recorded by the ordinary receiversdeployed on the sea floor alone for marine surveys or on land for landsurveys. The transmitter data may in some cases not be reliable, or maynot even be available due to failures of the monitoring system.

In one embodiment of the invention, the current inventive methodincludes the following two basic steps:

(a) detecting anomalous transmitter-pulse location(s) or windowsdirectly from the raw ordinary-receiver data in the time domain; and

(b) determining and correcting phase errors caused by anomaloustransmitter pulses using the continuity of phase and its derivative(with respect to source-receiver offset) with source-receiver offsetand/or time-domain correlations of the receiver waveforms, frequencyscaling behavior at near offset.

One embodiment of the present inventive method can be more fullysummarized in the following description, which refers to the flow chartof FIG. 7

At step 71, data from an ordinary survey receiver is examined foranomalies that could affect the phase of later data. This may beaccomplished by scanning receiver-gather data in the time domain usingzero-crossing-interval detection and/or cross-correlation techniques(both explained below) to locate locations and/or time (or offset)windows exhibiting transmitter anomalies. Alternatively, this may beaccomplished by identifying anomalous phase variations (versus offset)of the data in the frequency domain. At step 72, one such anomaly isselected for phase correcting. At step 73, the relative time shiftacross the selected anomaly window is determined, possibly by using thecross-correlation technique to applied to the data in the time domain.Alternatively, a phase shift can be determined by phase matching thedata in the frequency domain. In the latter case, the receiver data aretransformed from time domain to frequency domain by using a commoninitial time or phase (may be arbitrary) for all receivers in onetowline. In the frequency domain, the phase shift across an anomalywindow can be determined by matching the phase of the receiver data at afrequency ω at the bin just after the anomaly to the phase at the binjust before the anomaly. (Decomposition into the frequency domaindivides the data into components, each component associated with aparticular frequency in the source's frequency spectrum.) It may benoted that this phase error, Δφ, may be equated to a correspondingtiming error (the same timing error that could be determined directly inthe time domain) by the relationship Δφ=ωΔt. In practice, this phaseshift may be observable on data from some receivers but not on data fromreceivers whose signals have fallen below the noise floor. At step 74,the anomaly is fixed (i.e., a correction is made) by applying a phaseshift correction exp(−iωΔt) to all data (including data below the noisefloor) just after the anomaly window by using the time shift (or phaseshift) determined at step 73.

At step 75, it may be decided to fix another anomaly discovered at step71. Steps 72-74 are typically repeated until all transmitter anomaliesdetected in step 71 for the chosen receiver gather in the chosen towlineare corrected. When the time duration of a transmitter anomalousvariation is too long for the phase match technique to work in step 73,the cross-correlation technique can be used to determine the time shiftin the time domain. Typically, but not necessarily, when a time shift isdetermined in the time domain in step 73, the correction in step 74 isapplied also in the time domain. Similarly, when a phase shift isdetermined in the frequency domain in step 73, the correction in step 74is conveniently applied also in the frequency domain, one frequencycomponent at a time.

At this point in the present inventive method, the phase along atowline, or phase plotted against offset, has been made continuous forthe selected receiver gather and towline. But there still may be aglobal phase or time error (for example, an error in the choice of thetime origin, t=0), or there may be a global phase shift of the wholephase curve (for an individual frequency) along the tow line (phase vs.offset curve). A preferred way to check for whether a global correctionis needed, and to determine the appropriate correction (step 76), is touse the phase-scaling behavior at near offset, a technique disclosed byWillen (PCT International Patent Application No. PCT/US06/46329). Thistechnique is performed on the data components in the frequency domain.

In typical applications, steps 71-76 are repeated for all (or selected)receiver gathers on the chosen tow line, and for all (or selected) towlines in the survey (step 77).

The time-domain pulse measured by an ordinary receiver is a filteredversion of the transmitter pulse. The filtering effect comes from theelectromagnetic wave propagation in the earth's subsurface. Thevariation in pulse shape from one period to the next should be verysmall due to the small distance traveled by the transmitter in onetransmitter-pulse period, normally less than a few meters. Any abruptchanges of pulse shape, including the pulse period, are indications oftransmitter-signal anomalies. Those anomalies can be detected from theraw ordinary-receiver data as long as the signal is above the noisefloor.

Transmitter anomalies can be detected by zero-lag cross correlationbetween two consecutive bins of the receiver data in time domain. Byperforming cross correlation between two bins on both sides of atransmitter anomaly, the relative time shift across the anomaly can bedetermined. Once the time shift across an anomaly is determined in timedomain by the cross-correlation technique, the anomaly can be correctedby using a different initial transmitter phase corresponding to the timeshift to the data after the anomaly in the transformation from timedomain to frequency domain. This is equivalent to applying afrequency-dependent phase shift of exp(−iωΔt) across an anomaly infrequency domain. The time shift can also be corrected in time domain bydirectly shifting the time axis to the data after the anomaly.

The phase error caused by a transmitter anomaly can also be corrected infrequency domain by phase matching. Matching the data phase at one sideof a transmitter anomaly to that at another side is based on the factthat the phase of the CSEM data from a transmitter should be continuousversus the transmitter-receiver offset. In most practical applications,the first-order derivative of the phase with respect to offset can alsobe treated to be continuous versus offset. This can be used toextrapolate the phase across a time (or offset) window of anonymoustransmitter signals to correct the phase error introduced by thetransmitter anomalies.

The application of the phase correction to transmitter anomalies withina towline makes the phase consistent along the whole towline. The phasefrequency-scaling technique can fix the global time shift error of thetowline and any time shifts caused by un-detected and un-correctedtransmitter anomalies that occurred in time windows where the signalsfall below the noise floor by using near-offset data of the towline.When the time or phase errors are caused by transmitter anomalies, thedetected time shift across an anomaly from one receiver gather can beused to correct all receiver gathers sharing the same towline across thesame transmitter anomaly. The combination and repeat application ofthose techniques can correct most of the phase errors caused bytransmitter anomalies in a CSEM survey without using thesource-monitoring data.

Some receiver gathers may still have undetermined global time errorsafter the process when there are no near-offset data to correct theglobal time error of the corresponding towline. The global time errorscannot be determined accurately from the ordinary receiver data aloneunder this situation without further assumption on subsurface geology.But, these situations are rare. A more detailed description of theinvention follows.

Detecting Anomalous Variations of Transmitter Signal

As shown in FIG. 2, the pulses recorded by a channel at a receiverlocation have the same period and similar shape if the transmittersignal behaves normally as designed. An anomalous transmitter-signalvariation at a particular transmitter location, such as that at 21 inFIG. 2D, will cause an anomalous pulse on the receiver data, as shown inFIGS. 2A-C. Thus, one can scan the data recorded by ordinary receiversto determine locations (or times) of any anomalous pulses. This can bedone, for simple transmitter waveforms such as a square wave, by simplydetermining and plotting the time intervals between neighboringzero-crossings. For more complicated transmitter waveforms, thistechnique may fail due to the presence of too many zero crossings withina period or the absence of clear zero crossings. In this case, zero-lagcross correlation of two successive time bins (each having a width ofone or more periods of the received waveform), will show the presence ofanomalies at locations where the zero-lag cross-correlation coefficientsare reduced compared to locations of normal transmitter signals. If ananomalous zero-crossing time interval (or zero-lag cross-correlationcoefficient) appears on all of the three channels at the same time(which corresponds to the same transmitter location), this anomaly ismost likely caused by an anomalous transmitter pulse.

As stated, FIGS. 2A-D show a single-pulse anomaly, with a half periodlonger than the designed value of 4 s, recorded by each of the threechannels. This anomaly is detected directly from receiver data on thezero-crossing time-interval plots shown in FIGS. 4A-C by the singlespike 41 on the right side of the receiver location 42 in each of thethree channels. In this instance, the anomaly can be confirmed from thetransmitter monitoring data displayed in FIG. 2D. The transmitter cansometimes be unstable for a period of time much longer than one periodof the pulse, as detected and shown in FIGS. 4A-C by the multiple spikesat the left side of the receiver location 42. Each spike indicates atime interval between zero crossings different from 4 s, i.e., ananomalous pulse. The tails at both ends of each channel indicate thearea where the signal falls below the noise floor.

FIGS. 5A-D show the detection of anomalous transmitter waveforms by adifferent technique: the cross-correlation technique. FIG. 5B shows aplot of the zero-lag cross-correlation coefficient of two successivebins versus time (or bin) corresponding to the receiver data shown inFIG. 5A. The locations showing cross-correlation coefficients less thana critical value (the critical value being less than but close to 1.0)indicate the presence of abnormal transmitter pulses. The same schemecan also be applied to the recorded transmitter data, if available, forcross checking of the transmitter variations.

Anomalous source pulses can also be spotted after transforming thereceiver signal into the frequency domain. Rapid amplitude and phasevariations and phase discontinuities of the receiver data versustransmitter-receiver offset can be observed at offset windows where thetransmitter pulses are anomalous. In the frequency domain, the offsetrange affected by the anomalous variation of transmitter signal isusually larger than the actual duration of the anomalous transmittervariation due to the finite time window size used in the transformationfrom the time to frequency domain. FIG. 6A shows a transmitter-signalanomaly in the frequency domain as indicated by dashed lines witharrows. An anomalous transmitter pulse can be detected from the ordinaryreceiver data as long as it is above the noise floor in at least one ofthe receiver gathers in the survey. The cross correlation and frequencydomain techniques for identifying transmitter anomalies will next beexplained in more detail.

Determining Transmitter Anomalies and Tune Shifts by Cross Correlation

Cross correlation is a process that measures how much two time series ofnumbers resemble each other. The process is explained, for example, atpage 18 in Yilmaz's treatise, Seismic Data Processing, Society ofExploration Geophysicists (1987). The cross correlation function C(τ) oftime delay τ between two functions ƒ(t) and g(t) is defined by

${{C(\tau)} = {\sum\limits_{t}{{f(t)} \cdot {{g( {t + \tau} )}/\lbrack {\sqrt{\sum\limits_{t}{{{f(t)} \cdot {f(t)}}*}}\sqrt{\sum\limits_{t}{{{g(t)} \cdot {g(t)}}*}}} \rbrack}}}},$where the * represents the complex conjugate of a function. In thepresent case, ƒ(t) and g(t) are receiver signals in two time windows(bins) to be cross-correlated. Both ƒ(t) and g(t) are assumed to beperiodic in the calculation of the cross-correlation function C(τ) withperiod equal to the number of time samples in the corresponding windowor bin. When the time lag τ is set to zero, the cross correlation iscalled zero-lag correlation.

FIGS. 5A-D illustrate a procedure by which the cross-correlation ofadjacent time windows (bins) may be used to detect anomalies and computetime shifts. The horizontal labels on each drawing are time measured inbin number. The horizontal axis on FIG. 5C is a composite' scale. Theaxis labels individual cross-correlation functions as 1, 2 . . . 24. Foreach curve, an excursion to the right represents a positive correlationand an excursion to the left represents a negative correlation. Forexample, curves 5, 7, 8, 9, and 17 are nearly flat, as the correlationsthey represent are very small. By contrast; correlations such as 1through 4 range from low (roughly −1) to high (nearly +1) values,indicating strong correlation between waveforms. A waveform may compriseone or several periods of transmitter pulses. The cross-correlationscale of each correlation (−1 to +1) is not shown in FIG. 5C. For any ofthe 24 cross-correlation curves in FIG. 5C, the amplitude correspondingto a time shift of zero is the zero-lag cross-correlation coefficientbetween the two bins associated with the cross-correlation curve.

FIG. 5A shows synthetic receiver data in the time domain. The dataconsists of a series of waveforms of decreasing amplitude, withanomalies occurring between approximately t=5 and t=10 and between t=15and t=20. FIG. 5B plots the normalized zero-lag cross-correlationcoefficient of adjacent waveforms (i.e., between two successive bins);as expected, when the waveform shape does not change from one period tothe next, the cross correlation attains its maximum value of one.However, for the anomalous time windows, the cross correlation becomessmall, and this behavior can be used to identify transmitter-waveformanomalies.

The individual cross-correlation functions between the successivewaveforms and certain selected reference waveforms are shown in FIG. 5C.Examining these can clearly differentiate the locations of transmitterphase anomalies from stable locations. The cross-correlation functionsbetween each successive waveform and the reference waveform at t=1,which is the first waveform, attain peak values of one at zero time lagfor the first four periods, and then drop to small correlation levelsfor the next five periods, in the anomalous zone. The next sixcross-correlation functions between the corresponding waveforms and thereference waveform at t=4, which is the last stable waveform before theanomalous window, show a maximum correlation at a time shift of around−13 samples (see FIG. 5D, which plots the time shift of the correlationpeak). This indicates the time shift between the current stable waveformand the last stable waveform, and is the amount of time shift Δt thatmust be applied to the data after t=10. Similarly, the smallcross-correlation values between t=15 and t=20, from the crosscorrelation between the corresponding waveforms and the referencewaveform at t=15, which is again the last stable waveform before theanomaly, indicate a transmitter-anomaly window. The time lag of themaximum cross correlation (about −22 samples) shows the time shiftneeded with respect to the reference waveform at t=15. FIG. 5C hashorizontal broken lines at time shifts of 13 and 22 samples, to helpshow how the data for FIG. 5D between bins 10 and 15 and again betweenbins 20 and 24 can be obtained from FIG. 5C. A time sample, the unit forthe vertical axis in FIGS. 5C and 5D, represents a discretization of thetime scale wherein data collected in the time domain are sampled at aselected regular time interval. One bin contains many samples. Thus,careful examination of FIG. 5A shows that the curve is a succession ofclosely spaced, discrete data points, or sample values:

While both the zero-lag cross correlation and the lag of maximum crosscorrelation depart from stable values in the presence of transmitteranomalies, in this example the zero-lag cross correlation is a morereliable indicator of an anomaly. That is, the anomalies are morereadily identified in FIG. 5B than in FIG. 5D. For example, t=9 has asmall lag value (FIG. 5D), which could indicate a relatively stablewindow, but its zero-lag cross-correlation value is small (FIG. 5B). Infact, FIG. 5A shows that this area of the data has no signal pulse andis dominated by noise. Therefore, the zero-lag cross correlation is amore robust indicator of anomalous zones.

The reference waveform is defined above as the earlier in time of thetwo waveforms that are cross correlated. In FIGS. 5C and 5D, thereference waveform was chosen to be the last stable waveform before ananomaly, and it was kept the same until another anomaly is encountered.The way the reference waveforms are chosen in FIG. 5C is for the purposeof illustrating the time shift between the two stable time windowsbefore and after an anomaly window, and represents but one possibleembodiment of, the present inventive method. In field CSEM data, thewaveform varies with time (or transmitter location) from the earth'sfiltering effect. The reference waveform is preferably chosen from atime as close to the time of the waveform it correlates with as possiblein order to avoid the waveform variation form the earth's filteringeffect. Once the time shift across an anomaly is determined, thereference waveform should preferably be changed to the receiver waveformfrom the time bin immediately prior to the time bin being tested in timeintervals of stable transmitter signals. When an anomaly has beendetected, the waveform in the last stable time bin prior to theanomalous zone is preferably used as the reference waveform in order todetermine the relative time shift across the anomaly.

In one embodiment of the invention, the cross-correlation technique isused to detect and correct for transmitter anomalies as follows:

1. Scan through the received data for each ordinary receiver channel,computing the cross-correlation function for adjacent bins of the datacorresponding to one (or more) periods of the transmitter signal. Whilethe zero-lag cross-correlation remains above a pre-set critical value(close to one), continue to update the reference waveform to be thewaveform from the previous time bin.

2. When the normalized zero-lag cross-correlation coefficient becomesless than the critical value, designate this time (t*) as that of ananomaly. Store the earlier (in time) of the two data bins used in thecross-correlation computation. This bin is the last uncorrupted receiverdata segment occurring before the anomaly, and becomes the new referencewaveform.

3. Continue scanning through the receiver data, cross correlating thereference waveform stored in step (2) with following data segments. Whenthe maximum cross-correlation regains a value above the critical value,the end of the anomalous region has been found. The time lag at whichthe maximum cross-correlation value was found gives the time shift Δtthat must be applied to all of the data after the time of the anomaly(t*) found, in step (2). Apply this time shift to all data after t*.Replace the reference waveform with the data in the time bin just afterthe end of the anomalous region, which was found in this step.

4. Repeat steps (1)-(3) until all data from the current receiver gatherhave been processed.

5. All transmitter-caused phase anomalies have been identified andcorrected for this receiver. Anomalous zones may be discarded, asdesired, in further processing steps.

It may be noted elsewhere that steps 1-3 above correspond to steps 71-73of the flowchart in FIG. 7 for the embodiment of the invention in whichwaveform cross correlation is used.

Determining Transmitter Anomalies and Time Shifts by Phase Matching

In the frequency domain, both amplitudes and phases in a common-receivergather can be plotted versus the transmitter-receiver offset (or thetransmitter location relative to the receiver) for each of the frequencycomponents. FIG. 6A shows the phase variation versustransmitter-receiver offset of a common-receiver gather at two differentfrequencies f1 (curve 61) and f2 (curve 62). The phase and its slopewith respect to offset should be continuous and slowly-varying withoffset in the offset windows where the transmitter signal is stable (seeleft and right parts of each of the curves in FIG. 6A). However, thephase or its slope will be discontinuous across a transmitter-signalanomaly location or window, such as the anomaly that occurs in thevicinity of offset 63 at the middle part of the curves in FIG. 6A. Thephase continuity with offset can be used to determine the phase error(or time shift) caused by transmitter-signal anomalies.

To do this, the phase at one side of the anomaly is matched to the phaseat the other side by applying a phase shift equal to the phase jump 64or 65 as indicated in FIG. 6A to one side of the data. The phase becomescontinuous across the anomaly location after the shift as shown in FIG.6B. (The rapid phase oscillations near the anomaly location areartifacts from the transformation from time domain to frequency domain.)The phase shift Δφ frequency f1 is related to a time shift Δt by therelationship Δt=Δφ/(2πf1). The phases of all frequencies should bematched across the anomaly by a single time shift if the anomaly iscaused be the transmitter. The time shift derived from the phasematching at one frequency is applicable for all other frequencies. Thishas great value for applications on real data. The lowest-frequency datahave high signal-to-noise ratio over an extended offset range comparedto data at higher frequencies. The time shift can be determined reliablyfrom data at the lowest frequency and can then be applied to data at allother frequencies.

The phase matching procedure described above can be applied to atransmitter anomaly of duration much larger than the period of a normalpulse by using the continuity of both the phase and its derivative withrespect to offset. The phase at one side of the transmitter anomaly canbe extrapolated to the other side by using its slope with respect tooffset.

The relative time shift across a time (or offset) window with anomaloustransmitter-signal variation can also be obtained from the receiver datain time domain by using the cross-correlation technique used fordetecting the transmitter anomalies. This time shift can be applied tothe data at one side of the anomaly achieving the phase match before andafter the anomaly. This technique is more effective than thephase-matching technique discussed above for transmitter anomalies ofextended duration due to its independence on the phase continuityproperty. When using the cross-correlation method to determine the phaseerror, the time shift determination is also performed in the timedomain. The correction steps may be performed in either the time orfrequency domain. When the transmitter anomaly duration is long, thecross-correlation method may be preferred over the phase-matching methodas phase continuity can be difficult to apply across large timeintervals.

The techniques for determining the phase error described above can beapplied to every detected anomaly along a towline. After finishing allthe phase matches, the phases along the towline become consistent.

When a transmitter-signal anomaly can be observed on several receivers,the time shift determined from one receiver gather may be appropriatefor the same anomaly appearing on other receivers, because the clocktiming errors and errors associated with a receiver are either correctedalready or are usually small compared to errors from atransmitter-signal anomaly. This is useful when the time shift of atransmitter anomaly can be determined accurately from one receiver, butis hard to determine from other receivers because of low signalamplitudes measured by those receivers.

The time shift determined in step 76 of FIG. 7 for a given receiverincludes both timing errors of the clock and transmitter anomalies notdetected in step 71 in time windows where the data fall below the noisefloor. The time shifts determined in step 76 from different receivers ina towline should be very close to a common time shift value, whichcorresponds to the initial phase error of the towline, if alltransmitter-signal anomalies have been detected and corrected in steps71-75 and errors from other sources are negligible. If some transmitteranomalies have occurred, but have not been detected in step 71, betweenany two receivers in a towline, the time shifts determined in step 76from the two receiver gathers will be different. Applying the timeshifts derived from step 76 to the corresponding receiver gathers willcorrect any transmitter anomalies that might have happened but haveneither been detected in step 71 nor been corrected in step 74.

The foregoing application is directed to particular embodiments of thepresent invention for the purpose of illustrating it. It will beapparent, however, to one skilled in the art, that many modificationsand variations to the embodiments described herein are possible. Forexample, typically a CSEM survey is performed over water, and the sourceis towed along a line (called the tow line) over stationary receivers.But the invention applies to surveys over land as well, where the sourcewill not need to be towed, but will nevertheless typically be moved todifferent offsets along a line. The line of different offsets may bereferred to generally as a source line. A typical survey will consist ofmany source lines, often a grid of parallel and perpendicular lines. Allsuch modifications and variations are intended to be within the scope ofthe present invention, as defined in the appended claims.

The invention claimed is:
 1. A method for correcting source-relatedphase errors in electromagnetic data signals recorded by receivers in acontrolled-source electromagnetic survey of a subsurface region,comprising: (a) selecting a source line from the survey and an ordinarysurvey receiver associated with the selected source line, and obtainingthe electromagnetic data measured by the selected receiver for sourcepositions along the selected source line, said obtained data beingreferred to as a “receiver gather”; (b) finding a source-related anomalyin the receiver gather by examination of said data; (c) determining fromthe receiver gather a time shift or corresponding phase shift due to theanomaly; and (d) using the time/phase shift to adjust the phase of thedata in the receiver gather; wherein one or more of (a)-(d) areperformed using a computer.
 2. The method of claim 1, furthercomprising: (e) repeating steps (b)-(d) until any other source-relatedanomalies in the selected receiver gather are found and corrected byphase adjustment.
 3. The method of claim 2, further comprising: (f)repeating steps (a)-(e) for selected other ordinary receivers associatedwith the selected source line.
 4. The method of claim 3, furthercomprising repeating steps (a)-(f) for selected other survey sourcelines.
 5. The method of claim 2, further comprising determining whetherthe receiver gather has a global time/phase error, and if so,deteiiiiining the error and correcting the gather by applying atime/phase shift.
 6. The method of claim 1, wherein an anomaly is foundby zero-lag cross-correlation between consecutive bins of the data inthe receiver gather in time domain.
 7. The method of claim 1, wherein ananomaly is found by determining and comparing time intervals betweenneighboring zero crossings when the gather data is plotted vs. time. 8.The method of claim 1, wherein an anomaly is found by transforming thereceiver gather to components in the frequency domain and looking for adiscontinuity when data phase for a component is plotted vs.source-receiver offset.
 9. The method of claim 1, wherein the time shiftdue to an anomaly is determined in time domain by performing crosscorrelation between two bins of the data in the receiver gather oneither side of the anomaly.
 10. The method of claim 1, wherein the phaseshift due to an anomaly is determined in frequency domain by matchingdata phase at the anomaly's beginning to that at the anomaly's end. 11.The method of claim 10, wherein the first order derivative of phase withrespect to source-receiver offset is also matched across the anomaly.12. The method of claim 5, wherein the global error is determined infrequency domain using a frequency scaling property of data phase atnear offsets.
 13. The method of claim 1, wherein the data in thereceiver gather is adjusted for the anomaly by applying a phase shiftcorrection of exp(−iωΔt) in frequency domain to all data just after ananomaly window, phase shift Δφ and timing error Δt being related byΔφ=ωΔt.
 14. The method of claim 1, wherein gather data below noise levelis adjusted based on a time/phase shift determined from data in the samereceiver gather that are above noise level.
 15. The method of claim 1,wherein the source-related anomaly in the receiver gather is found byexamining the data gather for changes in pulse shape or pulse period tooabrupt to have been caused by geological conditions.
 16. A method forproducing hydrocarbons from a subterranean region, comprising: (a)obtaining electromagnetic field data from a controlled-sourceelectromagnetic survey of the subterranean region; (b) obtainingprocessing of the survey data, said processing having corrected the databy adjusting them for source-related phase errors by: (i) selecting asource line from the survey and an ordinary survey receiver associatedwith the selected source line, and obtaining the electromagnetic datameasured by the selected receiver for source positions along theselected source line, said obtained data being referred to as a“receiver gather”; (ii) finding a source-related anomaly in the receivergather by examination of said data; (iii) determining from the receivergather a time shift or corresponding phase shift due to the anomaly; and(iv) using the time/phase shift to adjust the phase of the data in thereceiver gather; wherein one or more of (i)-(iv) are performed using acomputer; and (c) producing hydrocarbons found in the subsurface regionat a location corresponding to a resistivity anomaly in a map of thephase-adjusted electromagnetic field data.